Supply & Demand of New Zealand Emission Units
This report investigates the supply and demand of New Zealand emission units, within New Zealand’s Emission’s trading scheme.
It investigates historical demand, historical supply and analysis the holdings of the current stockpile.
If you enjoy it, or find it interesting, please let the author by sending them a Linkedin message.
The data used to generate this report is free and publicly available.
You can access the data, and the notebook used to generate this report, from this repository.
Enjoy!
Historical Emissions
Figure 1 show’s New Zealand’s historic gross and net emissions.
Biogenic methane dominates the emission footprint making up ~49% of total emissions.
This data was sourced from Stats New Zealand.
The ETS Emissions is an approximation and was calculated as the difference between net and biogenic methane emissions.
This simple method suggests there are an annual, 20 million tonnes of CO2, captured by the ETS.
Code
# plot historic emissions
raw_emissions = pd.read_csv("./raw_data/historical_emissions.csv")
# Fix the biogenic methane issue
biogenic_methane = raw_emissions[["year", "biogenic_ch4"]].drop_duplicates()
# Step 2: Create a new DataFrame with 'Biogenic CH4' as the 'key'
biogenic_ch4_rows = biogenic_methane.assign(key='Biogenic CH4', emissions=biogenic_methane['biogenic_ch4']).drop(columns='biogenic_ch4')
# Adding back to historic emisisons
historic_emissions = pd.concat([raw_emissions, biogenic_ch4_rows], ignore_index=True)
historic_emissions.drop(columns=['biogenic_ch4'], inplace = True)
historic_emissions_pivot = historic_emissions.pivot(index='year', columns = 'key', values = 'emissions').reset_index()
historic_emissions_pivot['ETS_Emission_Gap'] = historic_emissions_pivot['GHG emissions: net'] - historic_emissions_pivot['Biogenic CH4']
def label_every_n_points(data, n=5):
# Create a new list of labels, with every nth label kept and the rest set to None
return [label if index % n == 0 else None for index, label in enumerate(data)]
n = 5 # for example, label every 5th point
plot_columns = ['GHG emissions: gross',"GHG emissions: LULUCF"]
# Apply the function to your label list
gross_ghg_labels = label_every_n_points([f'{x:,.2f}M' for x in historic_emissions_pivot['GHG emissions: gross']], n)
gross_forestry_labels = label_every_n_points([f'{x:,.2f}M' for x in historic_emissions_pivot["GHG emissions: LULUCF"]], n)
net_ghg_labels = label_every_n_points([f'{x:,.2f}M' for x in historic_emissions_pivot['GHG emissions: net']], n)
biogenic_ch4_labels = label_every_n_points([f'{x:,.2f}M' for x in historic_emissions_pivot['Biogenic CH4']], n)
ets_emission_gap_labels = label_every_n_points([f'{x:,.2f}M' for x in historic_emissions_pivot['ETS_Emission_Gap']], n)
historic_emisisons_fig = pex.bar(data_frame=historic_emissions_pivot, x = 'year', y = plot_columns)
historic_emisisons_fig.add_trace(
go.Scatter(x=historic_emissions_pivot['year'], y=historic_emissions_pivot['GHG emissions: gross'], name="Gross GHG", mode = 'text',
text= gross_ghg_labels,
showlegend=False,
textposition='top center',
)
)
historic_emisisons_fig.add_trace(
go.Scatter(x=historic_emissions_pivot['year'], y=historic_emissions_pivot['GHG emissions: LULUCF'], name="Forestry", mode = 'text',
text= gross_forestry_labels,
showlegend=False,
textposition='top center',
)
)
historic_emisisons_fig.add_trace(
go.Scatter(x=historic_emissions_pivot['year'], y=historic_emissions_pivot['GHG emissions: net'], name="Net GHG", mode = 'lines+text',
text= net_ghg_labels,
textposition='top center',
)
)
historic_emisisons_fig.add_trace(
go.Scatter(x=historic_emissions_pivot['year'], y=historic_emissions_pivot['Biogenic CH4'], name="Biogenic Methane", mode = 'lines+text',
text= biogenic_ch4_labels,
textposition='top center',
)
)
historic_emisisons_fig.add_trace(
go.Scatter(x=historic_emissions_pivot['year'], y=historic_emissions_pivot['ETS_Emission_Gap'], name="ETS Emissions", mode = 'lines+text',
text= ets_emission_gap_labels,
textposition='top center',
)
)
historic_emisisons_fig.update_layout(title = "New Zealand Historical Emission Profile by source and year - Data from Stats NZ",
yaxis_title = "tCO2e (million)"
)
historic_emisisons_fig.write_image("historic_emissions.png", width = 1920, height=1080, scale =2)
historic_emisisons_fig.show()Demand
To assess the demand, we look at the 2023 reported emission data published by the EPA.
Neglecting agricultural emisison sources, we can gain an estimate of the demand, by company and sector.
Figure 2 shows the reported emissions for New Zealand non-agricultural emissions.
Note that this is an approximate for NZU ETS demand. There is likely false data. For instance, Air New Zeland may or may not have reported their international aviation flights, which do not currently carry an ETS liability.
Code
demand_df = pd.read_csv(os.path.join(data_folder, 'participant_emission_report.csv')).rename(columns={"Reported emissions\n (tCO_{2} e)": "Reported Footprint tCO2e"})
# remove agiculture
demand_df = demand_df[demand_df['Sector']!= 'Agriculture']
# Replace commas and convert the '-' to NaN
demand_df["Reported Footprint tCO2e"] = demand_df["Reported Footprint tCO2e"].str.replace(',', '')
demand_df["Reported Footprint tCO2e"] = demand_df["Reported Footprint tCO2e"].replace('-', np.nan)
# Convert the column to numeric, coercing errors to NaN
demand_df["Reported Footprint tCO2e"] = pd.to_numeric(demand_df["Reported Footprint tCO2e"], errors='coerce')
demand_df = demand_df.dropna(subset=["Reported Footprint tCO2e"])
emisison_threshold = 500000
emission_filter = demand_df["Reported Footprint tCO2e"]>= emisison_threshold
demand_fig = pex.bar(data_frame=demand_df[emission_filter], x='Participant', y="Reported Footprint tCO2e", color='Sector')
# Change to log axis
demand_fig.update_layout(title = "New Zealands Top Non-Agriculatural Emitters (>500 kTCO2e) <br><sub>Data obtained from the EPA for 2023 reported emissions<sub>",
yaxis_type="linear")
demand_fig.show()
demand_fig.write_image("emitter_demand_top_emitters.png", width = 1920, height=1080, scale =2)Figure 3 shows the emissions by sector.
- Liquid Fossil Fuels - fuels used for transport.
- Stationary Energy - fossil fuels used for stationary processes. Most commonly industrial heat, like boilers and power plants.
- Industrial Processes - These are industries that generate CO2 from the production of goods not asociated with the direct combustion of fossil fuels, such as the aluminium smelter, production of glass or gold. The EPA define Industrial Processes as emissions from **“metal, mineral or chemical transformations”.
- Waste - emisisons from landfills.
Liquid fossil fuels is the biggest source of emissions (~53%), followed by stationary energy (~40%).
Code
summary = demand_df.groupby("Sector")["Reported Footprint tCO2e"].sum().reset_index()
fig = pex.pie(summary, values="Reported Footprint tCO2e", names="Sector", title = "2023 Reported Emisisons by Sector <br><sub>Data obtained by from the EPA</sub>")
fig.show()
fig.write_image("emission_source_pie.png", width = 1920, height=1080, scale =2)Top emitters, neglecting agriculture.
| Participant | Sector | Activity | Schedule | Year | Reported Footprint tCO2e |
|---|---|---|---|---|---|
| Loading... (need help?) |
Supply
NZUs are generated into the ETS by three mechanisms:
- Removals
- Industrial Allocations
- Auctions
Removals
Figure 4 shows the annual allocation per removal activity.
Forestry dominates the supply of removal units.
Mandatory reporting every 5 years, leading to the 5 year spikes of removal units awarded.
Without forestry, there is approximately 3M NZUs generated annually through other removal mechanisms as shown in Figure 5. This is important because all forestry units have a liability for the emission units, in the event the forest is destroyed or lumbered. That is, the property owner of the forest must surrender NZUs in the event they lumber the forest, or it’s destroyed in an adverse weather event, like a fire or cyclone.
This shows that there is approximately 3M liabililty free units generated from removal activity.
Code
# Removal Activity
removal_df = pd.read_csv(os.path.join(data_folder, "removal.csv"))
removal_df = removal_df.rename(columns={"The_Year":"Year", "Number_Of_Units":"NZUs"})
removal_df.fillna(value="Forestry", inplace=True)
removal_bar = pex.bar(data_frame=removal_df, x='Year', y="NZUs", color= "Activity")
removal_bar.update_layout(title_text = "NZUs awarded by removal activity by year",
annotations=[
dict(
xref='paper',
yref='paper',
x=0.21,
y=1.01,
xanchor='center',
yanchor='bottom',
text='Data from Emissions Trading Register',
showarrow=False,
font=dict(
size=12
)
)
]
)
aggregated_data = removal_df.groupby('Year')['NZUs'].sum().reset_index()
#line graph
removal_bar.add_trace(
go.Scatter(x=aggregated_data['Year'], y=aggregated_data['NZUs'], name="Total NZUs", mode = 'text',
text= [f'{x:,.2f}M' for x in aggregated_data['NZUs'] / 1e6],
textposition='top center',
showlegend=False
)
)
removal_bar.write_image("historic_removal.png", width = 1920, height=1080, scale =2)
removal_bar.show()Code
# Assuming removal_df is your DataFrame and it's already defined
non_forestry_removal = removal_df[removal_df["Activity"] != "Forestry"]
# Create the scatter plot
non_forestry_scatter = pex.bar(data_frame=non_forestry_removal, x="Year", y="NZUs", color="Activity", text_auto=True)
non_forestry_scatter.update_layout(title = "Removal activity without NZU liability by activity type and time")
# Show the plot
non_forestry_scatter.show()
# non_forestry_scatter.write_image("historic_removal_sans_forestry_bar.png", width = 1920, height=1080, scale =2)Industrial Allocation
Similar to the European Union Emissions Trading Program, there is an industrial allocation, designed to support trade exposed industries and mitigate emission leakage.
The calculation of industrial allocation under the Climate Change Response Act in New Zealand is calculated as:
Final Allocation Entitlement for the Year
\[FA = LA \times \sum(PDCT \times AB)\]
Where:
- $FA$ = Final allocation entitlement for the year.
- $LA$ = Level of assistance for the activity for the year.
- $PDCT$ = Amount of each prescribed product from the eligible industrial activity produced in the preceding year.
- $AB$ = Prescribed allocative baseline for the product.The level of assistance \(LA\) decreases annually after 2020 due to a phase-out rate.
Figure 6 shows the historical allocaiton of units via the industrial allocation.
Code
industrial_allocation_barplot = pex.bar(df_grouped, x="year", y="industrial allocation", text_auto=True, title = "Industrial allocatations by reporting year")
industrial_allocation_barplot.update_layout(yaxis_title = "NZUs")
industrial_allocation_barplot.show()
industrial_allocation_barplot.write_image("industrial_allocation.png", width = 1920, height=1080, scale =2)Auctions
NZUS are distributed through quarterly auctions. The Figure 7 shows the reduction in auction base supply.
This neglects the cost containment reserve volume.
Code
# Your data
years = ['2023', '2024', '2025', '2026', '2027', '2028']
supply = [15, 14.1, 12.6, 10.7, 9.1, 7]
percentage_change = [None, -6.00, -10.64, -15.08, -14.95, -23.08] # None for the first year as there's no change
# Create subplots and mention that we want to share the x-axis
fig = make_subplots(specs=[[{"secondary_y": True}]])
# Add bar graph for NZU Supply on primary y-axis
fig.add_trace(
go.Bar(x=years, y=supply, name='NZU Supply (MT)', text=supply),
secondary_y=False,
)
# Add line graph for Percentage Change on secondary y-axis
fig.add_trace(
go.Scatter(x=years, y=percentage_change, name='Percentage Change', mode='lines+markers+text', text=percentage_change),
secondary_y=True,
)
# Add titles and labels
fig.update_layout(
title_text="NZU Auction Supply excluding cost containment reserve by year"
)
# Set x-axis title
fig.update_xaxes(title_text="Year")
# Set y-axes titles
fig.update_yaxes(title_text="NZU Supply (Million)", secondary_y=False)
fig.update_yaxes(title_text="Year on Year Change (%)", secondary_y=True)
# Show the figure
fig.show()
fig.write_image("auction_volumes.png", width = 1920, height=1080, scale =2)Stockpile
There is a stockpile of New Zealand Emission Units.
Data was obtained from the Environmental Protection Agency.
Figure 8 shows the historic stockpile by NZU type. Figure 9 shows the historic stockpile neglecting all forestry units - the liability free stockpile.
- Other account holders are participants without NZU liabilities. They are speculators participating in the market.
- Participants are forestry owners and businesses that must submit emission reports.
Overall there is a surpless of 160M NZUs. Of which, ~100 m are forestry units and carry a liability to the land owner in the event the trees are lumbered or destroyed and not replanted.
If the supply of NZUs were stopped, this accounts for 8 to 3 years of NZU demand if foresty units are sold or not.
Code
stockpile_df = pd.read_csv(os.path.join(data_folder, 'stockpile_data.csv')).rename(columns={"Holdings of NZU_AUC":"AuctionUnits"})
def convert_cells_to_integers(df:pd.DataFrame, columns:list):
# Iterate over all columns in the DataFrame
for col in columns:
# Check if the column is of object type (string)
if df[col].dtype == 'object':
# Replace commas and convert the column to integers
df[col] = df[col].str.replace(',', '').astype(int)
return df
col_list = [col for col in stockpile_df.columns if col != 'Date']
df = convert_cells_to_integers(df=stockpile_df, columns=col_list)
df["ForestryNZUs"] =df['Holdings of NZU_FA'] + df['Holdings of NZU_FE and NZU_PFSI']
foo = df[['Date', 'All other NZUs', 'ForestryNZUs', 'AuctionUnits']]
fig = pex.bar(data_frame=foo, x='Date', y = ['All other NZUs', "ForestryNZUs", 'AuctionUnits'], text_auto=True)
# marker_trace = pex.scatter(data_frame=df, x='Date', y='Total holdings of NZUs', text_auto=True)
fig.add_trace(
go.Scatter(x=df['Date'], y=df['Total holdings of NZUs'], name="Total Holdings",
mode = 'text',
# text_auto=True,
text= [f'{x:,.2f}M' for x in df['Total holdings of NZUs']/1e6],
textposition='top center',
showlegend=False
)
)
# fig.update_traces(texttemplate='%{text::.2f}')
fig.update_layout(title = "30th June NZU balance by NZU Type",
yaxis_title = "NZUS")
fig.show()
# fig.write_image("Stockpile_with_forestry.png", width = 1920, height=1080, scale =2)Who is Holding What?
The Department of Environment (DoE) have information on the holdings by participant type.
I lumped all participants into two categories, those that hold emission liabilities like forestors or emitters (participants), and those who do not (other).
Figure 10 shows the historic holdings by NZU type and liability status.
Approximately 25% of units are held by parties that do not cary liabilities, such as speculators and financial institutions holding on behalf of emitters.
A minute detail on my reasoning
The DoE have a category of participants, they defined as “Recieved NZUS” who are not participants.
Their website defines these holders as:
“individuals and organisations that have received NZUs from the Crown, but who are not participants. This includes, for example, some owners of pre-1990 forests, participants in the Permanent Forest Sink Initiative (PFSI), and organisations receiving an Industrial Allocation”
I think this is silly because participants that receive industrial allocation need to report their emissions and surrender units. Similarily if a permanent forest sink initiative is destroyed, the owners are liable to pay the NZU debt or replant.
Code
# ACCOUNT HOLDERS
stacked_holdings_plot = pex.area(data_frame=filtered_data, x="Date", y="Holdings",
color = "NZU Type",
facet_col = 'Account holder category',
hover_data= ["NZU Type"],
title = "NZU Holdings by Type and Account Holder Category",
labels = {"Holdings": "NZU Holdings"},
template = "none")
stacked_holdings_plot.update_traces(opacity=0.7) # Adjust opacity
stacked_holdings_plot.write_image("NZU_holdings_forestry_assembled.png", width = 1920, height=1080, scale =2)
stacked_holdings_plot.show()